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The  three-dimensional  microstructure  of  an  SOFC  anode  is  quantified  using  a  dual  beam  focused  ion  beam 
scanning  electron  microscopy  (FIB-SEM)  system  equipped  with  an  energy  dispersive  X-ray  spectroscopy 
(EDX)  unit.  The  microstructure  of  the  Ni-YSZ  anode  is  virtually  reconstructed  in  a  computational  field 
using  a  series  of  acquired  two-dimensional  SEM  images.  The  three-phase  boundary  (TPB)  density  and 
tortuosity  factors  are  carefully  evaluated  by  applying  two  different  evaluation  methods  to  each  parameter. 
The  TPB  density  is  estimated  by  a  volume  expansion  method  and  a  centroid  method,  while  the  tortuosity 
factors  are  evaluated  by  a  random  walk  calculation  and  a  lattice  Boltzmann  method  (LBM).  Estimates  of 
each  parameter  obtained  by  the  two  methods  are  in  good  agreement  with  each  other,  thereby  validating 
the  reliability  of  the  analysis  methods  proposed  in  this  study. 

©  2009  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

A  solid  oxide  fuel  cell  (SOFC)  is  one  of  the  most  promising  energy 
conversion  devices  because  of  its  high  efficiency  and  fuel  flexibility 
[  1  ].  Long-term  durability  is  one  of  the  most  important  requirements 
for  the  practical  application  of  the  SOFC  system,  and  great  efforts 
have  been  made  to  develop  materials  with  high  stability  and  elec¬ 
trodes  with  optimal  microstructures.  Porous  Ni-zirconia  cermets 
such  as  Ni-yttria  stabilized  zirconia  (Ni-YSZ)  and  Ni-scandia  stabi¬ 
lized  zirconia  (Ni-ScSZ)  are  the  most  commonly  used  anode  materi¬ 
als  in  SOFCs,  because  their  thermal  expansion  coefficient  can  match 
that  of  the  electrolyte,  and  they  can  effectively  extend  reaction 
sites  (three-phase  boundary,  TPB)  with  high  electrocatalytic  activ¬ 
ity  [2,3].  It  is  widely  recognized  that  the  anode  microstructure  has 
a  significant  impact  on  cell  performance  as  well  as  cell  durability 
[4-8].  For  example,  effects  of  morphology  change  [7-9]  and  depen¬ 
dence  on  starting  materials  and  fabrication  methods  have  been 
reported  [10-12].  However,  the  quantitative  relationship  between 
the  anode  microstructure  and  the  polarization  resistance  is  not  fully 
understood.  Recently,  direct  measurements  of  3D  SOFC  electrode 
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microstructures  have  been  performed  by  focused  ion  beam  scan¬ 
ning  electron  microscopy  (FIB-SEM)  [13-17]  and  X-ray  computed 
tomography  (XCT)  [18].  Through  these  3D  measurements,  impor¬ 
tant  microstructural  parameters  such  as  TPB  length  and  tortuosity 
factors  can  be  obtained.  It  is  expected  that  these  new  methods  will 
provide  important  information  for  quantitatively  connecting  the 
electrode  microstructure  to  its  polarization  characteristics. 

It  seems,  however,  that  there  is  little  discussion  on  the  evalua¬ 
tion  methods  to  estimate  TPB  density  and  tortuosity  factors  using 
a  3D  dataset,  even  though  the  accuracy  of  the  estimation  depends 
on  the  evaluation  procedures.  For  example,  it  is  apparent  that  the 
total  TPB  length  would  be  overestimated  if  the  lengths  of  TPB  edge 
segment  are  simply  summed  because  of  the  inevitable  step-like 
pattern  of  voxel  edges.  Only  limited  descriptions  can  be  found  in 
literature  about  the  calculation  method  of  TPB  length  from  a  3D 
dataset.  Golbert  et  al.  [19]  counted  all  the  voxels  neighboring  a 
TPB  edge  and  divided  this  overall  number  by  four.  Suzue  et  al. 
[20]  counted  the  cubic  voxel  perimeter  and  assumed  that  TPB  was 
20%  smaller  than  the  total  value.  Smith  et  al.  [15]  and  Wilson  and 
Barnett  [21]  estimated  TPB  length  from  2D  images  by  stereology. 
The  advantage  of  using  a  3D  dataset  is  lost  in  this  method  and  the 
accuracy  of  the  TPB  length  estimated  by  stereology  should  also 
be  discussed.  Wilson  et  al.  [13]  and  Izzo  et  al.  [18]  evaluated  the 
tortuosity  factors  by  solving  Laplace  equation  with  finite  element 
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method  (FEM)  and  lattice  Boltzmann  method  (LBM),  respectively. 
Gostovic  et  al.  [14]  and  Smith  et  al.  [15]  calculated  tortuosity  by 
tracking  the  pore  center  locations  and  measuring  the  lengths  of  the 
paths.  As  the  values  of  TPB  length  and  tortuosity  factors  depend 
on  the  evaluation  procedures,  those  procedures  must  be  validated 
before  discussing  the  relationship  between  microstructural  param¬ 
eters  and  polarization  characteristics. 

In  the  present  study,  the  3D  microstructure  of  a  conventional 
Ni-YSZ  anode  is  quantified  by  dual  beam  FIB-SEM.  An  in-lens  sec¬ 
ondary  electron  detector  is  used,  which  provides  a  clear  contrast 
between  the  Ni  and  the  YSZ  phases  with  submicron  resolu¬ 
tion.  From  the  reconstructed  3D  structure,  various  microstructural 
parameters  such  as  volume  fraction,  TPB  length,  and  tortuosity 
factors  are  quantified.  We  pay  strong  attention  to  the  evaluation 
procedures  for  TPB  length  and  the  tortuosity  factors  to  ensure 
reliability  of  the  final  estimates  so  that  those  estimates  can  be 
used  to  correlate  microstructure  with  polarization  characteristics 
in  our  future  investigation.  Two  different  analyzing  methods  are 
compared  for  the  acquisition  of  each  parameter.  TPB  density  is  esti¬ 
mated  by  the  volume  expansion  method  and  the  centroid  method 
proposed  in  this  study  while  the  tortuosity  factors  are  validated  by 
the  random  walk  calculation  and  the  LBM-based  calculation. 

2.  Experimental 

2.1.  Sample  preparation 

In  this  study,  we  examined  the  Ni-YSZ  cermet  anode 
(Ni:YSZ  =  50:50  vol.%,  YSZ:8  mol%Y203-Zr02)  of  conventional  but¬ 


ton  cell,  Ni-YSZ|YSZ|LSM.  The  anode  material,  NiO-YSZ,  was  mixed 
with  polyethylene  glycol,  screen  printed  on  a  YSZ  electrolyte,  and 
sintered  at  1400  °C  for  5h.  The  (La0.8Sro.2)o.97Mn03  (abbreviated 
as  LSM)  cathode  was  also  mixed  with  polyethylene  glycol  to  form 
slurry.  It  was  then  screen  printed  on  the  other  face  of  the  elec¬ 
trolyte  and  sintered  at  1 150  °C  for  5  h.  These  button  cells  could  be 
used  in  power  generation  experiments  performed  in  laboratories. 
However,  no  power  generation  test  was  conducted  on  a  cell  exam¬ 
ined  in  this  study.  After  the  anode  of  the  test  cell  was  reduced  at 
1000°C,  the  cell  temperature  was  lowered  to  room  temperature 
and  the  anode  was  supplied  to  the  FIB-SEM  observation  as  a  sam¬ 
ple.  The  sample  was  infiltrated  with  epoxy  resin  (Marumoto  Struers 
KI<)  under  vacuum  conditions  so  that  the  pores  of  the  porous  elec¬ 
trode  could  be  easily  distinguished  during  SEM  observation.  Cured 
sample  was  polished  using  an  Ar-ion  beam  cross-section  polisher 
(JEOL  Ltd.,  SM-09010)  and  made  available  for  the  FIB-SEM  (Carl 
Zeiss,  NVision  40)  observation. 

2.2.  FIB-SEM  observation 

Observation  and  quantification  of  the  3D  microstructure  of  the 
Ni-YSZ  anode  are  facilitated  by  using  Kyoto  University’s  FIB-SEM, 
NVision  40,  equipped  with  a  Gemini  FE-SEM  column  (Carl  Zeiss),  a 
zeta  FIB  column  (SUNT),  and  a  multichannel  gas  injection  system 
(SUNT).  The  FIB-SEM  system  also  has  EDX  and  electron  backscat- 
tered  diffraction  (EBSD)  units. 

Fig.  1  schematically  shows  a  typical  setting  for  FIB-SEM  obser¬ 
vation.  In  this  system,  two  beams  have  a  coincident  angle  of  54°. 
In  this  study,  an  in-lens  secondary  electron  detector  with  a  typical 
acceleration  voltage  of  1  -2  kV  was  used  for  carrying  out  the  obser¬ 
vation.  A  sample  electrode  is  set  up,  as  shown  in  Fig.  1.  The  front 
part  of  the  target  volume  is  removed  by  FIB  milling  before  observa¬ 
tion.  By  using  the  gas  injection  system,  carbon  is  deposited  on  the 
surface  of  the  target  volume  to  protect  it  from  undesired  milling 
and  to  prevent  the  charging  up  of  the  observation  surface. 

Fig.  2(a)  shows  the  SEM  image  of  a  sample  electrode  after 
appropriate  sample  preparation.  FIB-SEM  observation  proceeds  as 
follows:  The  surface  of  the  observation  area  is  slightly  milled  in  the 
z  direction  by  FIB  so  that  a  new  x-y  observation  surface  is  exposed 
for  SEM  imaging.  By  automatically  repeating  this  ion  milling  and 
SEM  imaging  sequence,  known  as  “cut-and-see”  operation  in  the 
NVision  40  system,  a  series  of  SEM  images  necessary  for  3D  struc¬ 
ture  analysis  is  acquired.  Because  image  drifting  is  an  unavoidable 
problem  during  image  acquisition,  and  precise  image  alignment  is 
essential  for  the  quantitative  analysis  of  the  3D  microstructure,  one 
additional  procedure  is  required  before  beginning  the  milling  and 
imaging  sequence.  To  assist  the  latter  alignment  process  of  SEM 


Fig.  2.  (a)  SEM  image  of  sample  anode  before  “cut-and-see”  operation  sequence  and  (b)  reference  marks  created  on  the  carbon  coating. 


Table  1 

Sizes  and  spatial  resolutions  of  three  samples. 
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X 

y 

z 

Sample  1 

Sample  dimension  (p,m) 
Number  of  voxels 

26.095 

981 

10.906 

410 

4.74 

79 

Sample  2 

Sample  dimension  (p,m) 
Number  of  voxels 

25.722 

967 

11.624 

437 

6.572 

106 

Sample  3 

Sample  dimension  (p,m) 
Number  of  voxels 

26.341 

994 

10.768 

407 

6.048 

84 

images,  reference  marks  are  created  on  the  carbon  layer  deposited 
on  the  sample  surface  using  FIB,  as  shown  in  Fig.  2(b).  Sequential 
dataset  acquisition  by  FIB-SEM  is  carried  out  at  three  different  loca¬ 
tions  of  one  sample  anode.  The  sizes  and  spatial  resolutions  of  three 
datasets  are  summarized  in  Table  1. 


2.3.  3D  micro  structure  reconstruction 


Fig.  4.  Reconstructed  Ni-YSZ  anode  microstructure  (25.722  (Jim  x 
11.624  |jim  x  6.572  (Jim;  green:  Ni,  yellow:  YSZ). 


The  3D  microstructure  of  the  Ni-YSZ  anode  is  virtually  recon¬ 
structed  in  a  computational  field  using  2D  SEM  images  obtained 
via  FIB-SEM  observation.  SEM  images  are  first  aligned  based  on  the 
reference  marks  created  on  the  carbon  layer  deposited  on  the  sam¬ 
ple  surface  (see  Section  2.2).  We  found  that  controlling  FIB  milling 
precisely  to  maintain  the  milling  direction  perpendicular  to  the 
observation  surface  throughout  the  “cut-and-see”  sequence  in  this 
study  did  not  result  in  any  noticeable  image  inclination.  Therefore, 
image  alignment  was  carried  out  only  by  parallel  shifting  in  the  x 
and  y  directions. 

Separation  of  three  phases  is  carried  out  for  each  SEM  image. 
Fig.  2(b)  shows  that  the  pore  region  filled  with  epoxy  resin  is  eas¬ 
ily  distinguished  because  it  appears  as  a  thick  black  region.  On  the 
other  hand,  separation  of  the  solid  part  into  the  Ni  phase  and  the 
YSZ  phase  requires  careful  consideration.  EDX  analysis  is  performed 
at  the  beginning  or  the  end  of  FIB-SEM  observation  in  order  to 
correlate  brightness  with  two  solid  phases.  Fig.  3  shows  an  exam¬ 
ple  of  correlated  EDX  and  SEM  images.  The  Ni  and  YSZ  phases  can 
be  clearly  distinguished.  Correlation  obtained  is  applied  to  all  the 
images  of  the  same  data  series.  Separation  of  the  three  phases  is 
performed  semi-automatically,  with  some  manual  corrections  to 
finalize  the  process. 

After  the  alignment  and  the  phase  separation  processes,  the 
dataset  is  ready  for  3D  structure  reconstruction.  Fig.  4  shows  an 
example  of  Ni-YSZ  anode  microstructure  reconstructed  in  this 
study. 


3.  Quantification  results  and  discussion 

Geometric  parameters  of  the  sample  anode  are  evaluated  using 
3D  structure  data.  TPB  density  and  tortuosity  factors  are  important 
geometric  parameters  not  only  for  quantitatively  characterizing 
the  anode  performance  but  also  for  developing  reliable  models  for 
numerical  simulations.  It  seems,  however,  that  there  is  no  estab¬ 
lished  method  to  estimate  TPB  density  and  tortuosity  factors  for 
a  FIB-SEM  dataset,  even  though  the  accuracy  of  the  estimation 
depends  on  the  evaluation  procedures.  In  this  study,  we  apply  two 
different  evaluation  procedures  to  each  parameter  and  compare 
the  results  of  each  procedure  in  order  to  ensure  reliability  of  the 
final  estimate. 

3.1.  Volume  fractions 

The  3D  reconstruction  process  was  carried  out  for  three  sam¬ 
ple  datasets,  and  its  results  were  used  for  the  quantification 
study.  Table  2  summarizes  volume  fractions,  which  are  among  the 
most  fundamental  properties,  of  each  phase  of  the  three  samples. 
Because  the  sample  anode  is  fabricated  to  be  Ni:YSZ  =  50:50  vol.%, 
the  ratio  of  Ni  and  YSZ  volume  fractions  is  expected  to  be  unity. 
This  ratio  is  actually  calculated  to  be  1.09,  1.01,  and  0.94  for  the 
three  samples,  respectively.  We  found  that  this  variation  in  the 
values  of  the  abovementioned  ratio  was  not  negligible;  however, 


Fig.  3.  (a)  EDX  and  (b)  SEM  images  matching  for  separation  of  three  phases.  Blue  and  red  regions  correspond  to  Ni  and  YSZ  phases,  respectively,  in  the  EDX  image. 
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Table  2 

Volume  fractions  (%)  of  three  phases. 


Ni 

YSZ 

Pore 

Sample  1 

27.4 

25.1 

47.5 

Sample  2 

25.3 

25.1 

49.6 

Sample  3 

24.5 

26.0 

49.5 

it  was  still  within  the  acceptable  range.  The  fact  that  even  a  fun¬ 
damental  property  such  as  volume  fraction  suffers  such  variation 
implies  that  a  large  sample  size  is  preferred  for  the  quantitative 
analysis.  Increasing  the  sample  size,  however,  inevitably  low¬ 
ers  spatial  resolution.  The  balance  between  the  sample  size  and 
the  spatial  resolution  is  an  essential  problem  of  this  method, 
and  it  is  strongly  related  to  the  fabrication  process  of  cells.  In 
the  latter  sections,  we  consider  sample  2  to  be  a  representative 
sample.  It  was  resampled  to  have  a  62  nm  x  62  nm  x  62  nm  cubic 
voxel  structure,  because  some  of  the  following  analysis  meth¬ 
ods  require  such  a  voxel  structure.  The  resulting  sample  size  was 
18.600  [Jim  x  8.432  |xm  x  6.200  [xm. 

3.2.  Three-phase  boundary  length 

3.2.1.  Volume  expansion  method 

TPB  forms  lines  in  the  3D  reconstructed  field.  If  we  slightly 
expand  each  phase  outward  in  the  virtual  field,  the  overlapped 
regions  form  tube-like  volumes  that  contain  TPB  lines  inside.  Using 
this  method,  centerlines  of  those  tubes  are  taken  as  lines  that  rep¬ 
resent  TPBs,  and  their  lengths  are  measured.  It  is  worth  noting 
that  lines  obtained  through  this  method  theoretically  match  TPBs 
if  the  spatial  resolution  of  3D  reconstruction  is  sufficiently  high 
and  volume  expansion  is  limited  to  be  infinitely  small.  The  valid¬ 
ity  of  this  method  is  examined  by  applying  it  to  a  well-defined 
problem.  Fig.  5(a)  shows  a  geometry  adopted  for  a  test  calcula¬ 
tion.  We  assume  two  spheres,  each  representing  one  phase,  are 
placed  in  the  surrounding  third  phase.  When  the  two  spheres  are 
partially  overlapped  as  shown  in  the  figure,  a  circle-shape  TPB  is 
formed.  The  volume  expansion  method  is  applied  to  this  prob¬ 
lem  varying  its  space  resolution  and  the  obtained  TPB  lengths  are 
compared  with  the  theoretical  values.  Fig.  5(b)  shows  the  results 
of  the  test  calculation  corresponding  to  the  geometry  shown  in 
Fig.  5(a).  After  conducting  many  test  calculations  by  changing  rel¬ 
ative  positions  of  the  two  spheres,  we  concluded  that  this  method 
provides  reasonable  results  for  structures  that  have  a  characteristic 
length  greater  than  1 0  voxels.  TPB  density  estimated  by  the  volume 
expansion  method  is  2.487  [xm  [x m-3  for  the  sample  used  in  this 
study. 

3.2.2.  Centroid  method 

The  phase  of  each  voxel  is  assigned  as  either  Ni,  YSZ,  or  pore. 
If  the  neighboring  four  voxels  comprise  every  three  phases,  and 
the  phases  of  the  diagonal  voxels  are  not  the  same,  the  line  seg¬ 
ment  surrounded  by  the  four  voxels  is  defined  as  the  three-phase 
boundary.  Then,  triangles  are  defined  by  the  neighboring  three 
midpoints  of  the  three-phase  boundary  segments.  The  three-phase 
boundary  length  is  calculated  as  the  distance  between  the  cen¬ 
troids  of  these  triangles.  The  total  TPB  length  and  the  active  TPB 
lengths  in  the  x,  y,  and  z  directions  are  listed  in  Table  3.  The  total 
TPB  length  obtained  by  the  centroid  method  is  2.556  |xm  |xm-3. 
The  difference  between  the  volume  expansion  method  is  less  than 
3%. 

Although  present  Ni-YSZ  composition  (Ni:YSZ  =  50:50  vol.%)  is 
not  the  same  as  those  in  Refs.  [13,17],  total  TPB  densities  obtained 
by  the  present  methods  are  smaller  than  those  reported  in  Ref.  [13], 
4.28  [xm  [xm-3,  and  in  Ref.  [17],  10.58  and  9.36  |xm  [xm-3.  Present 
results  show  better  agreement  with  the  values  obtained  from  stere- 


Low-* - Space  resolution - ►  High 


Fig.  5.  Validation  of  TPB  length  obtained  by  applying  the  volume  expansion  method 
to  a  well-defined  problem:  (a)  test  geometry  set  for  the  validation  and  (b)  compar¬ 
ison  with  the  theoretical  value  at  various  space  resolutions. 

ology  [21],  i.e.,  2.73  and  2.50  |xm  [xm-3  for  Ni:YSZ  =  44:56  vol.%  and 
55:45  vol.%,  respectively. 

3.3.  Tortuosity  factor 

3.3.1.  Random  walk  calculation 

Tortuosity  factor  can  be  statistically  calculated  from  the  ran¬ 
dom  walk  process  of  nonsorbing  particles.  As  the  first  step  of  this 
method,  a  number  of  random  walkers  are  randomly  distributed  to 
the  pore  voxels.  Each  walker  randomly  chooses  one  of  the  neigh¬ 
boring  voxels  as  its  possible  location  in  the  next  time  step.  If  the 
selected  neighboring  voxel  represents  the  pore  part,  the  walker 
migrates  to  that  voxel.  If  the  selected  voxel  is  a  solid  voxel,  the 
walker  stays  at  the  current  voxel  and  waits  for  next  time  step.  In  this 
procedure,  neither  absorption  nor  desorption  is  taken  into  account. 
When  repeating  this  procedure,  the  mean  square  displacement  of 
random  walkers  is  calculated  as  follows: 

n 

<i r2(t )>  =  i^({xKO-Xi(0)}2  +  {y;(t)-yK0)}2 

1=1 

+  {Zi(t)-Zi(0)}2)  (i) 

where  n  is  the  total  number  of  random  walkers  and  ()  denotes  an 
ensemble  average.  Because  mean  square  displacement  (r2)  is  pro¬ 
portional  to  time,  the  diffusion  coefficient,  D,  of  nonsorbing  species 
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Table  3 

Three-phase  boundary  length  by  centroid  method. 
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Total  TPB 

Active  TPB  (x  =  0  gm:  electrolyte,  x=  18.6  gm:  current  collector) 
Active  TPB  (x  =  0  gm:  current  collector,  x  =  18.6  gm:  electrolyte) 
Active  TPB  (y  =  0  gm:  electrolyte,  y  =  8.43  gm:  current  collector) 
Active  TPB  (y  =  0  gm:  current  collector,  y  =  8.43  gm:  electrolyte) 
Active  TPB  (z  =  0  gm:  electrolyte,  z  =  6.2  gm:  current  collector) 
Active  TPB  (z  =  0  gm:  current  collector,  z=  6.2  gm:  electrolyte) 


TPB  length  (gm  (Jim-3) 

2.556(100%) 

1.539(60.2%) 

1.067  (41.7%) 
1.400(54.8%) 
1.399(54.7%) 
1.647(64.4%) 
1.659(64.9%) 


is  related  to  the  time-derivative  of  (r2),  as  follows: 

1  d(r2(t)) 

6  dt 


(2) 


The  mean  square  displacement  in  a  porous  medium,  (r2)pore* 
takes  a  lower  value  than  that  obtained  in  a  free  space  because  of  the 
obstruction  effects  of  solids.  Considering  that  the  effective  diffusion 
coefficient,  Deff,  is  expressed  as 


neff  \r  1  d(r2(t))Pore 
°  =VW6 - dt - 


(3) 


the  degree  of  reduction  is  measured  quantitatively  by  the  tortuosity 
factor  as  follows  [22,23]: 


_  D  d(r2(t))  ; d(r2(f))Pore 

Pore  Deff/VPore  dt  /  dt 


(4) 


Note  that  tortuosity  factor  is  different  from  tortuosity,  which  is 
defined  as  the  ratio  of  the  average  winded  pore  length  to  the  thick¬ 
ness  of  the  porous  material  [24].  Tortuosity  factor  is  much  more 
important  than  tortuosity  because  the  former  is  directly  related  to 
effective  diffusivity  or  conductivity. 

When  the  porous  medium  has  an  anisotropic  pore  structure, 
the  mean  square  displacement  (r2)  may  be  divided  into  directional 
mean  square  displacements,  (x2),  (y2),  and  (z2),  as  follows: 


n 

(X2(t)>  =  -  Xi(°))2  (5) 

1  =  1 
n 

<y2(t)>  =  ^(yi(t)-yi(0)}2  (6) 

i=l 

n 

<z2(t)>  =  ^(z,(t)-zj(0)}2  (7) 

1=1 


Because  Eqs.  (5),  (6)  and  (7)  are  obtained  in  three  directions,  tortu¬ 
osity  factors  have  three  different  values,  corresponding  to  the  x,  y, 
and  z  directions. 

If  there  are  isolated  pores,  pores  connected  to  the  computational 
domain  surface  boundaries  are  extracted  from  stacked  images  as 
“effective  pores”.  Isolated  pores  are  excluded  from  the  calculation 
because  they  do  not  contribute  to  the  gas  diffusion  and  may  cause 
an  overestimation  of  the  tortuosity  factor.  In  this  case  the  volume 
fraction  of  the  effective  pores,  Vpore  eff,  is  used  in  Eq.  (3)  and  the 
tortuosity  factor  is  calculated  as  follows: 

D  lye  d(r2(t))  /d(r2(t))P  ore  (on 

Pore  “  De ff/VPore  “  Vpore>eff  dt  7  dt 

In  order  to  obtain  an  accurate  value  of  the  tortuosity  factor,  a 
long  time  step  and  a  large  number  walker  are  necessary.  How¬ 
ever,  as  the  time  step  proceeds,  random  walkers  may  go  out  of 
the  stacked  images.  A  mirror  symmetric  boundary  condition  is 
applied  to  solve  this  problem,  because  it  guarantees  connectivity  of 
pore  structure  across  boundaries.  The  calculation  is  carried  out  for 
1 0,000,000  time  steps  with  1 00,000  random  walkers.  Data  obtained 


during  the  first  5,000,000  steps  are  omitted  in  the  calculation  of  the 
diffusion  coefficient  because  walkers  have  not  yet  experienced  the 
porous  solid  structure  at  the  early  stage  of  the  random  walk.  We 
use  “SIMD-oriented  Fast  Mersenne  Twister  (SFMT)”  as  a  high  per¬ 
formance  pseudorandom  number  generator  in  this  method  to  meet 
the  requirement  for  a  large-scale  random  number. 

The  accuracy  of  this  method  is  examined  by  applying  it  to  a 
well-defined  problem  shown  in  Fig.  6(a).  An  inclined  straight  path 
is  formed  in  the  computational  domain.  The  tortuosity  factor  in 
x-direction,  rx,  can  theoretically  be  calculated  in  the  ideal  case  as 
rx  =  1  +  tan2  0  where  0  is  the  inclination  angle.  rx  is  calculated  vary¬ 
ing  the  path  width,  W,  and  the  inclination  angle.  The  results  are 
compared  with  the  theoretical  values  in  Fig.  6(b).  The  discrepancy 
from  the  theoretical  value  becomes  large  at  a  large  inclination  angle 


0  [degrees] 


Fig.  6.  Validation  of  tortuosity  factors  obtained  through  the  random  walk  calcula¬ 
tion  for  a  well-defined  problem:  (a)  the  problem  settings  for  the  validation  and  (b) 
comparison  with  the  theoretical  value  under  various  conditions. 
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when  W  is  small.  We  can  see  a  tendency  that  the  tortuosity  fac¬ 
tor  is  overestimated  when  the  space  resolution  is  insufficient.  The 
requirement  for  the  space  resolution  is  severer  than  what  required 
for  the  TPB  length  estimation.  We  concluded  that  this  method  pro¬ 
vides  reasonable  results  for  structures  that  have  a  characteristic 
length  greater  than  20  voxels. 

The  calculation  results  of  this  method  are  shown  in  the  next 
section  together  with  the  results  of  LBM-based  calculation. 


3.3.2.  Lattice  Boltzmann  method  (LBM)-based  calculation 

Assuming  that  Ni  and  YSZ  are  perfect  electronic  and  ionic  con¬ 
ductors,  the  gaseous,  electronic,  and  ionic  diffusion  equations  are 
solved  inside  each  of  the  obtained  3D  structures  of  Ni,  YSZ,  and  pore 
phases  as  follows: 


_d_ 

dXa 

_d_ 

dxa 

d 

dXrv 


dc 


dxa 

^el  dcf e- 
F  dxa 

gjo 

2  F  dxa 


=  0 


=  0 


=  0 


(9) 


(10) 


(11) 


where  xa  represents  x,  y,  z  directions;  C  is  concentration;  D  is  dif¬ 
fusion  coefficient;  r]  is  electrochemical  potential;  o  is  conductivity; 
and  F  is  Faraday’s  constant.  The  lattice  Boltzmann  method  (LBM) 
[20]  is  used  to  solve  Eqs.  (9),  (10)  and  (11).  The  LB  equation  with 
the  LBGK  model  in  the  collision  term  is  expressed  as  follows: 


(12) 


/f(x  +  Cf  At,  t  +  At)  =fi(x,  t )  -  p[fi(x,  t)  -_/fq(x,  t)]  +  w,  A t 

In  Eq.  (12),/;  represents  the  velocity  distribution  function  with 
velocity  c2-  in  the  ith  direction,  and  ffq  is  the  Maxwellian  local 
equilibrium  distribution.  For  the  3D  LBM  simulation,  the  D3Q15 
(i  =  1  -1 5)  or  D3Q1 9  (i  =  1  -1 9)  models  are  commonly  used.  However, 
it  has  been  shown  that  in  case  of  simple  diffusion  simulation,  the 
D3Q6  (i  =  1-6)  model  can  be  used  with  only  a  slight  loss  of  accuracy 
[25];  therefore,  we  have  used  the  D3Q6  model  in  our  current  study. 
The  relaxation  time,  t *  =  0.99,  is  fixed  for  all  simulations.  Dirichlet 
boundary  conditions  are  applied  at  the  boundary  surfaces.  From  the 
LBM  calculation,  we  can  obtain  the  effective  diffusion  coefficient 
Deff  of  the  gas  phase  as  well  as  the  respective  effective  conductiv¬ 
ities  of  the  Ni  and  YSZ  phases,  cr^f  and  a^z.  Because  Eqs.  (9),  (10) 
and  ( 1 1 )  are  similar,  tortuosity  factors  for  each  phase,  rpore >  tNi  *  and 
rYsz>  defined  as  follows: 

neff 


eff 

Tel 


Teff 


 Wore 

(13) 

^Pore 

£ 

^1- 

II 

(14) 

m 

Wsz_ 

—  _  ^io 

(15) 

?YSZ 

where  Wore,  Wii,  and  Vysz  are  the  volume  fractions  of  the  Pore,  Ni, 
and  YSZ  phases,  respectively. 

Calculated  tortuosity  factors  for  the  three  phases  are  summa¬ 
rized  in  Table  4  together  with  the  values  estimated  by  the  random 
walk  method.  The  tortuosity  factors  evaluated  in  these  methods 
only  include  geometrical  aspects.  Other  effects  such  as  reaction  or 
diffusion  on  the  surface  of  electrochemically  active  nickel  are  not 
considered.  The  results  afforded  by  both  methods  are  in  good  agree¬ 
ment  with  each  other  (less  than  3%  difference).  The  values  for  the 
pore  in  Table  4  are  also  close  to  the  values  reported  in  Wilson  et 
al.  [13]  (rx,  Ty,  tz)  =  (2.1,  2.2,  1.9).  To  the  authors’  knowledge,  this 
is  the  first  report  on  the  tortuosity  factors  of  the  solid  phases,  rNj 
and  rYsz-  Cross-sectional  areas  normal  to  the  x,  y,  and  z  directions 


Table  4 

Anisotropic  tortuosity  factors  calculated  by  two  methods. 


Random  walk 

Lattice  Boltzmann  method 

rpore 

X 

2.05 

2.03 

y 

1.99 

2.06 

z 

1.78 

1.83 

rNi 

X 

22.10 

21.68 

y 

29.46 

29.45 

z 

6.91 

6.94 

*\sz 

X 

27.89 

27.66 

y 

14.95 

14.82 

z 

9.86 

9.84 

are  52.3  f±m2,  115.3  f±m2,  and  156.8  f±m2,  respectively.  As  shown 
in  Table  4,  the  tortuosity  factors  of  the  solid  phases,  rNi  and  tYSz, 
have  large  values  for  the  x  and  y  directions.  On  the  other  hand,  rpore 
has  nearly  the  same  values  for  all  the  three  directions.  We  conclude 
that  the  winded  electronic  and  ionic  paths  are  disconnected  at  the 
side  boundaries  for  the  x  and  y  directions.  From  Table  4,  it  is  evident 
that  the  processed  volume  size  is  not  sufficiently  large  for  evaluat¬ 
ing  the  effective  conductivities  of  the  solid  phases  for  the  present 
sample. 

The  validated  TPB  length  and  tortuosity  factors  will  be  coupled 
with  the  numerical  simulation  in  the  future  work,  which  will  pro¬ 
vide  a  quantitative  relationship  between  microstructure  and  anode 
polarization. 

4.  Conclusions 

The  three-dimensional  microstructure  of  a  Ni-YSZ  anode  was 
quantified  by  dual  beam  FIB-SEM.  An  in-lens  secondary  electron 
detector  was  used,  which  provided  clear  contrast  between  the  Ni 
and  YSZ  phases  as  well  as  submicron  resolution.  The  reconstructed 
three-dimensional  structure  was  used  to  quantify  microstructural 
parameters  such  as  the  volume  fraction,  TPB  length,  and  tortuos¬ 
ity  factors.  TPB  density  was  estimated  by  the  volume  expansion 
method  and  the  centroid  method  proposed  herein;  the  results 
matched  well— less  than  3%  difference  was  found  between  the  two 
sets  of  results.  This  proves  the  reliability  of  the  two  analysis  meth¬ 
ods  for  TPB  density  proposed.  The  tortuosity  factors  were  evaluated 
by  the  random  walk  calculation  and  the  LBM-based  calculation. 
The  results  for  the  three  phases  afforded  by  both  methods  matched 
well— differences  were  less  than  3%.  This  validates  the  two  analysis 
methods  for  the  tortuosity  factors  adopted  in  this  study.  However, 
it  was  also  revealed  that  the  volume  of  the  sample  used  in  this  study 
was  not  sufficient  for  a  quantification  of  the  effective  conductivi¬ 
ties  of  Ni  and  YSZ.  The  next  step  will  be  coupling  the  validated  TPB 
length  and  tortuosity  factors  with  the  numerical  simulation,  which 
will  provide  a  quantitative  relationship  between  microstructure 
and  anode  polarization. 
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